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ABSTRACT 

The forms of electron density structures in kinetic Alfven wave turbulence 
are studied in connection with scintillation. The focus is on small scales 
L ~ 10^ — 10^° cm where the Kinetic Alfven wave (KAW) regime is active in 
the interstellar medium, principally within turbulent HII regions around bright 
stars. MHD turbulence converts to a KAW cascade, starting at 10 times the ion 
gyroradius and continuing to smaller scales. T hese scales are inferred to d o minate 



scintillation in the theor y of Boldyrev et al. (iBoldyrev &: Gwinnl l2003al Jbl. 12005 



Boldyrev fc Koniglll2006l ). From numerical solutions of a decaying kinetic Alfven 
wave turbulence model, structure morphology reveals two types of localized struc- 
tures, filaments and sheets, and shows that they arise in different regimes of re- 
sistive and diffusive damping. Minimal resistive damping yields localized current 
filaments that form out of Gaussian-distributed initial conditions. When resis- 
tive damping is large relative to diffusive damping, sheet-like structures form. 
In the filamentary regime, each filament is associated with a non-localized mag- 
netic and density structure, circularly symmetric in cross section. Density and 
magnetic fields have Gaussian statistics (as inferred from Gaussian-valued kur- 
tosis) while density gradients are strongly non-Gaussian, more so than current. 
This enhancement of non-Gaussian statistics in a derivative field is expected since 
gradient operations enhance small-scale fluctuations. The enhancement of density 
gradient kurtosis over current kurtosis is not obvious, yet it suggests that modest 
fluctuation levels in electron density may yield large scintillation events during 
pulsar signal propagation in the interstellar medium. In the sheet regime the same 
statistical observations hold, despite the absence of localized filamentary struc- 
tures. Probability density functions are constructed from statistical ensembles in 
both regimes, showing clear formation of long, highly non-Gaussian tails. 



Subject headings: ISM: electron density — ISM: general — MHD — Turbulence 
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Introduction 



Models of scintillation have a long history. Many (ILee &: Jokipiilll975al Jbl: ISuttonlll97l[ ) 
carry an implicit or explicit assumption of Gaussian statistics, applying to either the electron 
density field itself or its autocorrelation function (herein referred to as "Gaussian Models"). 
Lee fc Jokipiil (Il975al ) is a representative approach. The statistics of the two-point correlation 
function of the index of refraction e(r), A{p) = J dz' {t{x,z)t{x + p,z')) determines, among 
other effects, the scaling of pulsar signal width r with dispersion measure DM. The index 
of refraction e(r) is a function of electron density fluctuation n(r). The quantity A{0) enters 
the equations, representing the second moment of the index of refraction. If the distribution 
function of e(r) has no second-order moment (as in a Levy distribution) A{0) is undefined. 
The assumption of Gaussian statistics leads to a scaling of r ~ DM^, which contradicts 
observation for pulsars with DM > 30 cm~'^ pc. For these distant pulsars, r ~ DM"^ 



flSuttonlll97ll : lBoldvrev fc Gwinnll2003al Jbf) 



To explain the anomalous DM^ scaling, ISuttonI (Il97l[ ) argued that the pulsar signal 
encounters strongly scattering turbulent regions for longer lines of sight, essentially arguing 
that the stati stics, as sarnpled by a pulsar signal, are nonstationary. Considering the pulse 
shape in time, IWilliamsonI ( 1l972l . Il973l Il974l ) is unable to match observations with a Gaussian 
Model of scintillation unless the scattering region is confined to 1/4 of the line of sight between 
the pulsar and Earth. These assumptions may have physical basis, since the ISM may not be 
statistically stationary, be ing composed of different regions with varying turbulence intensity 
teoldvrev fc GwinnI [20051 ). 



The theory of iBoldvrev fc GwinnI fl2003al Jbl. |2005|); iBoldvrev fc Konid (120061 ) takes a 
different approach to explain the anomalous DM^ scaling by considering Levy statistics for 
the density difference (defined below). Levy distributions are characterized by long tails with 
no defined moments greater than first-order [i.e., ^(0) is undefined for a Levy distribution]. 
The theory recovers the r ~ DM'^ relation with a statistically stationary electron density 
field. This theory also does not constrain the scattering region to a fraction of the line-of-sight 
distance. 

The determinant quantity in the theory of Boldyrev et al. is the density difference. 
An = n(xi, z) — n(x2, z). According to this model, if the distribution function of An has a 
power-law decay as |An| — »• oo and has no sec ond moment, then it is possible to recover the 
r ~ DM'^ scaling (IBoldyrev fc GwinnI l2003b( ). Assuming sufficiently smooth fluctuations. 
An can be expressed in terms of the density gradient, (t{z)\ n(xi) — n(x2) ~ a{z) ■ (xi — X2). 
Perhaps more direc t ly, the density gradient enters the ray tracing equations [Eqns. (7) in 
Boldyrev fc GwinnI ( l2003al )]. and is seen to be central to determining the resultant pulsar 
signal shape and width. This formulation of a scintillation theory does not require the 
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distribution of An to be Gaussian or to have a second-order moment. 



The notion that the density difference is characterized by a Levy distribution is a con- 
straint on dynamical models for electron density fluctuations in the ISM. Consequently the 
question of how a Levy distribution can arise in electron density fluctuations assumes con- 
siderable importance in understanding the ISM. 

Previous work has laid the groundwork for answering this question. It has been es- 
tablished that electron density fluctuations associated with interstellar magnetic turbulence 
unde rgo a significant ch ange in character near the scale lOp^, where ps is the ion sound gyrora- 
dius ( iTerry et al.ll200ll ). At larger scales, electron density is passively a dvected by the turbu- 
lent fl ow of an MHD cascade mediated by nonlinear shear Alfven waves (iGoldreich &: Sridhar 
19951 ). At smaller scales the electron density becomes compressive a nd the turbulent e nergy 
is carried into a cascade mediated by kinetic Alfven waves (KAW) (ITerry et al.ll200ll ). The 
KAW cascade brings electron density into equipartition with the magnetic field, allowing 
for a significant increase in amplitude. The conversion t o a KAW cascade has been ob- 
served in numerical solutions of the gyrokinetic equ ations (IHowes et al.ll2006l). and is consis- 



tent with observations from solar wind turbulence (IHarmon fc Coles 



2005 



Bale et al. 2005 



Leamon et al.lll998[ ). Importantly, it puts large amplitude electron density fiuctuations (and 



large amplitude density gradients) a t the gyroradius scale (^ 10 — 10 cm), a desirable set 



of conditions for pulsar scintillation (IBoldyrev fc Koniglll2006l ). It is therefore appropriate to 
consider whether large-amplitude non-Gaussian structures can arise in KAW turbulence. 

This question has been partially answered in a stu dy of current filament formation in 
decaying KAW turbulence (ITerry fc Smith! l2007l . |2008| ) . In numerical solutions to a two- 
field model with broadband Gaussian initial conditions large amplitude current filaments 
spontaneously arose. Each filament was associated with a large-amplitude electron density 
structure, circular in cross-section, that persisted in time. These electron density structures 
were not as localized as the corresponding current filaments, but were coherent and not 
mixed by surrounding turbulence. The observation of large amplitude current filaments is 
simila r to the large-ampl itude vortex filaments found in decaying 2D hydrodynamic turbu- 
lence ( IMcWilliamsl Il984j ). Counterparts of such structure s in 3D are predicte d to be the 
dominant component for higher order structure functions (IShe &: Levequelll994j ). 



Terry fc SmithI (120071 ) proposed that a nonlinear refractive magnetic shear mechanism 



prevents the structures from mixing with turbulence. Radial shear in the azimuthally directed 
magnetic field associated with each large-amplitude current filament decreases the radial 
correlation length of the turbulent eddies and enhances the decorrelation rate. Eddies are 
unable to persist long enough to penetrate the shear boundary layer and disrupt the structure 
core. The structure persistence mechanism allows large- amplitude fiuctuations to persist for 
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many eddy-turnover times. As the turbulent decays these structures eventually dominate the 
statistics of the system. The spatial structure of the density field associate d with localized 



circu larly symmetric current filaments was shown from analytical theory (ITerry &: Smith 



20071 ) to yield a Levy-distributed density gradient field. The kurtosis for the current field 
was significantly larger than the Gaussian- valued kurtosis of 3, indicating enhanced tails. 
The electron density and magnetic field kurtosis values were not significantly greater than 3. 
However, just as the current is non-Gaussian when the magnetic field is not, it is expected 
that numerical solutions should show non-Gaussian behavior for the density gradient. In the 
present paper, density gradient statistics are measured and found to be non-Gaussian. Rather 
than relying on kurtosis values alone, the probability density functions (PDFs) are computed 
from ensembles of numerical solutions, showing non-Gaussian PDFs for the density gradient 
field. 

The previous studies of filament generation in KAW turbulence leave significant unan- 
swered questions relating to structure morphology and its effect on scintillation. It is well 
established that MHD turbulence admits structures that are both filament-like and sheet- 
like. Can sheet-like structures arise in KAW turbulence? If so, what are the conditions or 
parameters favoring one type of structure versus the other? If sheet-like structures dominate 
in some circumstances, what are the statistics of the density gradient? Can they be suffi- 
ciently non-Gaussian to be compatible with pulsar scintillation scaling? It is desirable to 
consider such questions prior to calculation of rf wave scattering properties in the density 
gradient fields obtained from numerical solutions. 

In this paper, we show that both current filaments and current sheet structures naturally 
arise in numerical solutions of a decaying KAW turbulence model. Each has a structure of 
the same type and at the same location in the electron density gradient. These structures 
become prevalent as the numerical solutions progress in time, and each is associated with 
highly non-Gaussian PDFs. Moreover, we show that small-scale current filaments and current 
sheets, along with their associated density structures, are highly sensitive to the magnitude 
of resistive damping and diffusive damping of density fluctuations. Current filaments persist 
provided that resistivity r] is small; similarly, electron density fluctuations and gradients 
are diminished by large diffusive damping in the electron continuity equation. The latter 
results from collisions assuming density fluctuations are subject to a Pick's law for diffusion. 
The magnitude of the resistivity affects (1) whether current filaments can become large in 
amplitude, (2) their spatial scale, and (3) the preponderance of these filaments as compared 
to sheets. The magnitude of the diffusive damping parameter, /i, similarly influences the 
amplitude of density gradients and, to a lesser degree, influences the extent to which electron 
density structures are non-localized. 
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In the ISM resistive and diffusive damping become important near resistive scales. How- 



ever, it is well kn own that collisionless damping effects are also present (iLysak fc Lotkdll996 



Bale et al.ll2005[ ). and quite possibly dominate over collisional damping in larger scales near 



the ion Larmor radius. The collisional damping in the present work is understood as a 
heuristic approach that facilitates analysis of the effects of different damping regimes on 
the statistics of electron density fluctuations. By varying the ratio of resistive and diffusive 
damping we can, as suggested above, control the type of structure present in the turbulence. 
This allows us to isolate and study the statistics associated with each type of structure. It 
also allows us to assess and examine the type of environment conducive to formation of the 
structure. We consider regimes with large and small damping parameters, enabling us to 
explore damping effects on structure formation across a range from inertial to dissipative. 
Future work will address collisionless damping in greater detail. 



1.1. Background Considerations for Structure Formation 



The coherent structures observed in numerical solutions of decaying KAW turbulence, 
whether elongated sheet s or localized filaments, are sim ilar to structures observed in decaying 
MHD turbulence, as in iKinney &: McWilliamd (Il995[ ). In that work, the flow field initially 
gives rise to sheet-like structures. After selective decay of the velocity field energy, the system 
evolves into a state with sheets and filaments. During the merger of like-signed filaments, 
large-amplitude sheets arise, limited to the region between the merging filaments. These 
short-lived sheets exist in addition to the long-lived sheets not associated with the merger 
of filaments. In the two-field KAW system, however, there is no flow; the sheet and filament 
generation is du e to a different mech anism, of which the filament generation has previously 
been discussed (ITerry fc SmithI 120071 ). 



Other work (IBiskamp fc Welterl Il989l : iPolitano et al.l Il989[ ) observed the spontaneous 
generation of current sheets and filaments in numerical solutions, with both Orszag-Tang 
vortex and randomized initial conditions. These 2D reduced MHD numerical solutions mod- 
eled the evolution of magnetic flux and vorticity with collisional dissipation coefficients t], 
the resistivity, and u, the kinematic viscosity. The magnetic Prandtl number, z//?7, was set 
to unity. These systems are incompressible and not suitable for modeling the KAW system 
we consider here - they do however illustrate the ubiquity of current sheets and filaments, 
and serve as points of comparison. For Orszag- Tang-like initial conditions with large-scale 
flux tubes smooth in profile, current sheets are preferred at the interfaces between tubes. 
Tearing instabilities can give rise to filamentary current structures that persist for long times, 
but the large-scale and smoothness of flux tubes do not give rise to strong current filaments 
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localized at the center of the tubes. To see this, consider a given flux tube, and model it as 
cylindrically symmetric and monotonically decreasing in r with characteristic radial extent 



for < r < a. The current is localized at the center with magnitude 



J =-4. (2) 

Thus flux tubes with large radial extent a have a corresponding small current filament 
at their center. Hence, initial conditions dominated by a few large-scale flux tubes are 
not expected to have large amplitude current filaments at the flux tube centers, but favor 
current sheet formation and filaments associated with tearing instabilities in those sheet 
regions. At X points cu rrent sheet folding and filamentary structures can arise (see, e.g. 



Biskamp fc Welterlll989l . Fig. 10), but these regions are small in area compared to the quies- 
cent fiux-tube regions. Note that if, instead of Orszag- Tang-like initial conditions, the initial 
state is random, one expects some regions with flux tubes that have a small, and therefore 
a sizable current filament at the center. 

Consider now the effect of comparatively large or small 1]. In the case of large t], the 
central region of a flux tube is smoothed by the collisional damping, thus having a strong 
suppressive effect on the amplitude of the current filament associated with such a flux tube. 
Large-amplitude current structures are localized to the interfaces between flux tubes. In 
the process of mergers between like-signed filaments (and repulsion between unlike-signed) , 
large current sheets are generated at these inte rfaces, similar to the large-am plitude sheets 



generated in MHD turbulence during mergers (IKinney fc McWilliamsl Il995[ ). For small r], 
relatively little suppression of isolated current filaments should occur; if these filaments are 
spatially separated owning to the buffer provided by their associated flux tube, they can 
be expected to survive a long time and only be disrupted upon the merger with another 
large-scale flux tube. Large rj, then, allows current sheets to form at the boundaries between 
flux tubes while suppressing the spatially-separated current filaments at flux tube centers. 
Small r) allows interface sheets and spatially separated filaments to exist. 

These simple arguments suggest that the evolution of the large-amphtude structures and 
their interaction with turbulence is thus strongly influenced by the damping parameters. As 
such, the magnitudes of the damping parameters are expected to affect the resultant pulsar 
scintillation scalings. The present paper considers the effect of variations of these damping 
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parameters, rj and fi. In the KAW model, the (unnormahzed) resistivity takes the form 
1] = mf.Vf./ne^ and the density diffusion coefficient is = plvei where irtf. is the electron mass, 
Vf. is the electron collision frequency, n is the electron density, e is the electron charge, and 
Pe = VTel'^ce^ with VTe the electrou thermal velocity, and Uce the electron gyrofrequency. The 
ratio of these terms, c^rj/Airii = 2//3, where /3 = SnnkT / B'^ is the ratio of plasma to magnetic 
pressures. When we vary this ratio, as we will do in the numerical solutions presented here, 
we have in mind that we are representing regions of different /3. However, as a practical 
matter in the numerical solutions, we must vary the damping parameters independently of 
the variation of /3, since the kinetic Alfven wave dynamics require a small /3 to propagate. For 
the warm ioniz ed medium, typ ical parameters are Tg = 8000 K,n = 0.08 cm~'^, \B\ = 1.4 fiG, 



6B = 5.0 nG (jFerrierell200ll ). With these parameters, the plasma /3 formally ranges from 



0.05 — 1.2, spanning a range of plasma magnetization. 

We present the results of numerical solutions of decaying KAW turbulence to ascertain 
the effect of different damping regimes on the statistics of the fields of interest, in particular 
the electron density and elec tron density gradient. In the 77^/1 regime (using normalized 



parameters), previous work (jCraddock et al.l Il99ll ) had large-amplitude current filaments 



that were strongly localized with no discernible electron density structures {fj, was large to 
preserve numerical stability). This regime is unable to preserve density structures or density 
gradients. The numerical solutions presented here have 17 ~ /i and 1] ^ fi; in each limit 
the damping parameters are minimized so as to allow structure formation to occur, and are 
large enough to ensure numerical stability for the duration of each numerical solution. We 
investigate the statistics of both filaments and sheets in the context of scintillation in the 
warm ionized medium. 

The paper is arranged as follows: section |2] gives an overview of the KAW model and 
normalizations, its regime of validity, and its dispersion relation. Section [3] discusses the 
numerical method used and the field initializations. The negligible effect of initial cross- 
correlation between fields is discussed. Results for the two damping regimes are given in 
section |H where the type of structures that form, whether sheets-and-filaments or predomi- 
nantly sheets, are seen to be dependent on the values of t] and /i. PDFs from ensemble nu- 
merical solutions are presented in section illustrating the strongly non-Gaussian statistics 
in the electron density gradient field for both the 1] fi and r] ^ fi regimes. This suggests 
that non-Gaussian electron-density gradients are robust to variation in 77, as long as the 
overall damping in the continuity equation is not too large. Some discussions regarding the 
limitations of numerical approximation for this work and possible enhancements-particularly 
a model that addresses driven KAW turbulence-are given in concluding remarks. 
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2. Kinetic Alfven Wave Model 



The kinetic Alfven wave (KAW) model u sed in this paper is the sani e model used in the- 
ories of pulsar scintilla tion through the ISM ( iTerry &: Smith! 120071 . |2008| ) and in earlier work 
( ICraddock et al.lll99ll ). It is a reduced, two-field, small-scale limit of a more general reduced 
three-field MHD system (jHazeltinelll993l : iRahman Sz Wielandlll983l : [Fernandez fc Terrylll997l ) 
that accounts for electron dynamics parallel to the magnetic field. 

The 3- field model applies to large and small-scale fluctuations as compared to ps, the ion 
gyroradius evaluated at the electron temperature. In large-scale strong turbulence magnetic 
and kinetic fluctuations are in equipartition, with electron density passively advected. In 
the limit of small spatial scales (< lOps) the roles of kinetic and internal fluctuations are 
reversed - magnetic fluctuations are in equipartition with density fluctuations, and kinetic 
energy experiences a go-it-alone cascade without participating in the magnetic-internal en- 
ergy interaction. The shear- Alfven phy sics at large scale is supplanted by kinetic-Alfven 
physics at small scale (ITerry et al.ll200ll ). 



In the Boldyrev et al. theory, the length scales that dominate scintillation for pulsars 
with DM > 30 pc cm~^ are small, around 10* — 10^° cm. This motivates our focus on the 
small-scale regime of the more general 3-fleld system. The dominant interactions are between 
magnetic and internal fluctuations, via kinetic Alfven waves. In these waves, electron density 
gradients along the magnetic fleld act on an inductive electric fleld in Ohm's law. The electron 
continuity equation serves to close the system. The normalized equations are 



9,7/- = Vp + r/oJ-r^sW, (3) 

dtn = -V\\J + poV^n- fi2^^^^n, (4) 

V|| = (9, + W X z ■ V, (5) 

J = V^^, (6) 

with ifj = {Cs/c)eAz/Te, the normalized parallel component of the vector potential and 
n = {Cs/VA)n/no the normalized electron density. The normalized resistivity is r]o = 
{c^/4:hVaPs)Vsp, with rjsp the Spitzer resistivity, given in the introduction. The normal- 
ized diffusivity is p = plve/ PsVa- The time and space normalizations are ta = Ps/Va and 
Ps = Cs/^i- Here Cg = {Te/rriiY/'^ is the ion acoustic velocity, Va = B / {AuminQY/'^ is 
the Alfven speed, and f2j = eB/rriiC is the ion gyrofrequency. Electron density diffusion is 
presumed to follow Pick's law; more detailed damping would necessarily consider kinetic ef- 
fects and cyclotron resonances. The 7]2 and p2 terms (hyper-resistivity and hyper-diffusivity) 
are introduced to mitigate large-scale Fourier-mode damping by the linear diffusive terms. 



- 9 - 



Throughout the remainder of the paper, we drop the 2 subscript from ri2 and fj,2 and refer 
to the hyper-dissipative terms as r] and fi. 

Three ideal invariants exist: total energy E = J d'^x[{'Vip)'^ + an'^]; flux F = j d^x if)"^ 
and cross-correlation He = J d'^x nip. Energy cascades to small scale (lar ge k) while the flux 



and c ross-correlation undergo an inverse cascade to large-scale (small k) (IFernandez fc Terry 



19971 ). The inverse cascades require the initialized spectrum to peak at fco 7^ to allow for 



buildup of magnetic flux at large-scales for later times. 

Linearizing the system yields a (dimensional) dispersion relation u = VAkzk±ps- The 
mode combines perpendicular oscillation associated with a finite gyroradius with fluctuations 
along a mean field (z-direction). The oscillating quantities are magnetic field and density, 
out of phase by 7r/2 radians. 

In the limit of strong mean field, quantities along the mean field (2;-direction) equilibrate 
quickly, which allows d/dz — )■ 0, or k^ — 0. Kinetic Alfven waves still propagate, as long 
as there are a broad range of scales that are excited, as in fully developed turbulence. As 
kz — )■ 0, all gradients are localized to the plane perpendicular to the mean field. Presuming a 
large-scale fluctuation at characteristic wavenumber ko, smaller-scale fluctuations propagate 
linearly along this larger-scale fluctuation so long as their characteristic scale k satisfies 
k ^ kf). In this reduced, two dimensional system, the above dispersion relation is modified 
to be w = VA(bfcQ ■ k/B)kps which is still Alfvenic but with respect to a perturbed large-scale 
amplitude perpendicular to the mean field. Relaxing the scale separation criterion yields 
uj (X k"^ for the general case. 



Numerical Solution Method 



We evolve Eqs. ([3]) and (jl]) in a 2D periodic box, size [27r] x [27r] on a mesh of resolution 
512 X 512. The ip and n scalar fields are evolved in the Fourier domain, with th e nonlinearitie s 
advanced pseudospectrally and with full 2/3 dealiasing in each dimension (jOrszagl Il97ll ) . 
The diffusive and resistive terms normally introduce stiffness into the equations; using an 
integrating factor remo ves any stabil i ty con straints stemming from these terms. Following 
the scheme outlined in ICanuto et al.l ( 1l990l ). we start with the semi-discrete formulation of 
Eqs. (ED and dH): 



dt 



(7) 
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duk 
dt 



-fik^nk-J^[V\\J] 



where J^[-] denotes the discrete Fourier transform. We do not exphcitly expand the nonhnear 
terms as they will be integrated separately. The hyper-damping terms (proportional to k^) 
are included above. Damping terms corresponding to the Laplacian operator (proportional 
to /c^) are not included in this section for clarity, but are trivial to incorporate. Equations [7] 
and [S] can be put in the form 



d 
~dt 



d r 
'dt 



-e"'' [V||J] . 

A second-order Runge-Kutta scheme for the V'fc difference equation is 



(9) 
(10) 



m+l/2 



-r,k^At/2 \,m 



(11) 



with a similar form for the scheme. 



m+l/2 



(12) 



3.1. Initial conditions 

The 'j/'A: and fields are initialized such that the energy spectra are broad-band with a 
peak near /cq ~ 6 — 10 and a power l aw spectrum for k > k n. The falloff in k is predicted 
to be for small-scale turbulence. ICraddock et al.l ( 119911 ) use k~^, between the current- 
sheet limit of k~^ and the kinetic- Alfven wave strong-turbulence limit of The numerical 
solutions considered here have either k~^ or k^^. The only qualitative difference between the 
two spectra is the scale at which structures initially form. The k'"^ spectra has more energy 
at smaller scales, leading to smaller characteristic structure size. After a few tens of Alfven 
times these smaller-scale structures merge and the system resembles the initial k~^ spectra. 

The Hk and ipk phases can be either cross-correlated or uncorrelated. By cross-correlated 
we mean that the phase angle for each Fourier component of the Uk and ipk fields are equal. 
In general. 
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n^=\A,,\e'\ ijk = \Bk\e''\ (13) 

where \Ak\ and \Bk\ are the Fourier component's amphtude, set according to the spectrum 
power-law. For cross-correlated initial conditions, 6i = 62 for all k at the initial time. 
For uncorrelated initialization, there is no phase relation between corresponding Fourier 
components of the Uk and ipk fields. 



Craddock et al.l (119911 ) focused on the formation and longevity of current filaments in 
a turbulent KAW system. To preserve small-scale structure in the current filaments, these 
numerical solutions set 77 = and had /i ~ 10~^, with a resolution of 128 x 128, corresponding 
to a kmax of 44. Large-amplitude density structures that would have arisen were damped 
to preserve numerical stability up to an advective instability time of a few hundred Alfven 
times, for the parameter values therein. The numerical solutions presented here explore 
a range of parameter values for rj and fi. They make use of hyper-diffusivity and hyper- 
resistivity of appropriate strengths to preserve structures in n, B and J. An advective 
instability is excited after ~ 10^ Alfven times if resistive damping is negligible. The 77 = 
solutions-not presented here due to their poor resolution of small-scale structures-see large- 
amplitude current filaments arise, but they can be poorly resolved at this grid spacing. With 
no resistivity, the finite number of Fourier modes cannot resolve arbitrarily small structures 
without Gibbs phenomena resulting and distorting the current field. 

We have found through experience that small hyper-resistivity and small hyper-diffusivity 
preserve large-amplitude density structures and their spatial correlation with the magnetic 
and current structures, while preventing the distortion resulting from poorly-resolved current 
sheets and filaments. They allow the numerical solutions to run for arbitrarily long times, 
and the effects of structure mergers become apparent. These occur on a longer timescale 
than the slowest eddy turnover times. The results presented here will consider two regimes 
of parameter values, the rj ^ ^ and rj ^ fi regimes. The effect of cross-correlated and 
uncorrelated initial conditions will be addressed presently. 



4. Results 

It is of interest to examine whether cross-correlated or uncorrelated initial conditions 
affect the long-term behavior of the system. Two representative numerical solutions are 
presented here that reveal the system's tendency to form spatially-correlated structures in 
electron density and current regardless of initial phase correlations. This study establishes 
the robustness of density structure formation in KAW turbulence and lends confidence that 
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such structures should exist in the ISM under varying circumstances. The first numerical 
solution has cross-correlated initial conditions between the n and ip fields; the second, un- 
correlated. Damping parameters rj and /i are equal and large enough to ensure numerical 
stability while preserving structures in density, current and magnetic fields. These examples 
also serve to explore the intermediate 77/// regime. 

The energy vs. time history for both numerical solutions are given in Figs. [1] and |2l 
Total energy is a monotonically decreasing function of time. The magnetic and internal 
energies remain in overall equipartition throughout the numerical solutions. Magnetic en- 
ergy increases at the expense of internal energy and vice versa. This energy interchange is 
consistent with KAW dynamics and overall energy conservation in the absence of resistive 
or diffusive terms. The exchange is crucial in routinely producing large amplitude density 
fluctuations in this two-fleld model of nonlinearly interacting KAWs. 

The total energy decay rates for the uncorrelated and correlated initial conditions in 
Figs. [1] and [2] differ, with the latter decaying more strongly than the former. The damping 
parameters are identical for the two numerical solutions, and the decay-rate difference re- 
mains under varying randomization seeds. The magnitudes of the nonlinear terms during 
the span of a numerical solution in Eqs. (|3]) and dl]) for uncorrelated initial conditions are 
consistently larger than those of correlated initial conditions by a factor of 5. This difference 
lasts until 2500 Alfven times, after which the decay rates are roughly equal in magnitude. 
The steeper energy decay during the run of numerical solutions with uncorrelated initial con- 
ditions (Fig. 12]) suggests that the enhancement of the uncorrelated nonlinearities transports 
energy to higher k (smaller scale) more readily than the nonlinearities in the correlated case. 
Relatively more energy in higher k enhances the energy decay rate as the linear damping 
terms dissipate more energy from the system. The initial conflguration, whether correlated 
or uncorrelated, is seen to have an effect on the long-term energy evolution for these de- 
caying numerical solutions. It will be shown below, however, that the correlation does not 
signiflcantly affect the statistics of the resulting flelds. 

For cross-correlated initial conditions, we expect there to be a strong spatial relation 
between current, magnetic fleld and density structures through time. Figs. [3] and S] show 
the n and |B| contours at various times. For the latest time contour, the spatial structure 
alignment is evident. Further, in Fig. [5l the circular magnetic fleld structures (magnetic fleld 
direction and intensity indicated by arrow overlays) align with the large-amplitude density 
fluctuations. The correlation is evident once one notices that every positive-valued circular 
n structure corresponds to counterclockwise-oriented magnetic fleld, and vice versa. Fig. O 
is at a normalized time of 5000 Alfven times, deflned in terms of the large Bq. The system 
preserves the spatial structure correlation indeflnitely, even after structure mergers. 



- 13 - 



The second representative numerical solution is one with uncorrelated initial conditions. 
Contour plots of density and |B| are given in Figs. [6] and [71 respectively. It is noteworthy that, 
similar to the cross-correlated initial conditions, spatially correlated density and magnetic 
field structures are discernible at the latest time contour. 

In Fig. [8] the circular density structures correspond to circular magnetic structures. 
Unlike Fig. [5] the positive density structures may correspond to clockwise or counterclock- 
wise directed magnetic field structures. This serves to illustrate that, although the initial 
conditions have no phase relation between fields, after many Alfven times circular density 
structures spatially correlate with magnetic field structures and persist for later times. 

The kurtosis excess as a function of time, defined as -ft'(S) = (H^) /o"| — 3, is shown 
in Figs. M and [TU] for correlated and uncorrelated initial conditions, respectively. Positive 
K indicates a greater fraction of the distribution is in the tails as compared to a best-fit 
Gaussian. These figures indicate that the non-Gaussian statistics for the fields of interest 
are independent of initial correlation in the fields. In particular, the density gradients, |Vn|, 
are significantly non- Gaussian as compared to the current. Because scintillation is tied to 
density gradients, this situation is expected to favor the scaling inferred from pulsar signals. 

The tendency of density structures to align with magnetic field structures regardless of 
initial conditions indicates that the initial conditions are representative of fully-developed 
turbulence. After a small number of Alfven times the memory of the initial state is removed 
as the KAW interaction sets up a consistent phase relation betwee n the fiuctuations in the 



magnetic and density fields. Previous work ( iTerry fc Smith! 120071 ) presented a mechanism 



whereby these spatially correlated structures can be preserved via shear in the periphery of 
the structures. The above figures indicate that this mechanism is at play even in cases where 
the initial phase relations are uncorrelated. 

In the damping regime presented above, circularly symmetric structures in density, 
current and magnetic fields readily form and persist for many Alfven times, until disrupted 
by mergers with other structures of similar amplitude. It is possible to define, for each 
circular structure, an effective separatrix that distinguishes it from surrounding turbulence 
and large-amplitude "sheets" that exist between structures, [see, e.g., the magnetic field 
contours at later times in Fig. [71] The density field has significant gradients in both the 
regions surrounding the structure and within the structures themselves. The ability to 
separate these circular structures from the background sheets and turbulence is determined 
by the magnitudes - relative and absolute - of the damping parameters. Larger damping 
values erode the small-spatial-scale structures to a greater extent and, if large enough, disrupt 
the structure persistence mechanism that, for a fixed diameter, depends on a sufficiently large 
amplitude current filament to generate a sufficiently large radially sheared magnetic field. 
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The preceding results were for a damping regime where t]/^ ~ 1, an intermediate 
regime. Numerical solutions with yU = and t] small explore the regime where i]/^ — t- 0. 
In this regime, which is opposite the regime used in Craddock et al., circularly symmetric 
current and magnetic structures are not as prevalent, rather, sheet-like structures dominate 
the large amplitude fluctuations. Current and magnetic field gradients are strongly damped, 
and the characteristic length scales in these fields are larger. 

Contours of density for a numerical solution with /x = are shown in Fig. [TTJ Visual 
comparison with contours for runs with smaller damping parameters (Fig. El where rj = ^) 
indicate a preponderance of sheets in the /i = case, at the expense of circularly-symmetric 
structures as seen above. All damping is in t]; any current filament that would otherwise form 
is unable to preserve its small-scale, large amplitude characteristics before being resistively 
damped. Inspection of the current and \B\ contours for the same numerical solution [Figs. 
[T2] and [13] reveal broader profiles and relatively few circular current and magnetic field 
structures with a well-defined separatrix as in the small r] case. Since there is no diffusive 
damping, gradients in electron density are able to persist, and electron density structures 
generally follow the same structures in the current and magnetic fields. 

Kurtosis excess measurements for the n = numerical solutions yield mean values 
consistent with the = /i numerical solutions, as seen in Fig. [TH Magnetic field strength and 
electron density statistics are predominantly Gaussian, with current statistics and density 
gradient statistics each non-Gaussian. Perhaps not as remarkable in this case, the density 
gradient kurtosis excess is again seen to be greater than the current kurtosis excess - this 
is anticipated since the dominant damping of density gradients is turned off. W i th few er 



filamentary current structures, however, the mechanism proposed in lTerry fc Smith! (120071 ) is 
not likely to be at play in this case, since few large-amplitude filamentary current structures 
exist. Sheets, evident in the density gradients in Fig. [T5] and in the current in Fig. [Tulare the 
dominant large-amplitude structures and determine the extent to which the density gradients 
have non-Gaussian statistics. The current and density sheets are well correlated spatially. 
The largest sheets can extend through the entire domain, and evolve on a longer timescale 
than the turbulence. Sheets exist at the interface between large-scale flux tubes, and are 
regions of large magnetic shear, giving rise to reconnection events. With rj relatively large, 
the sheets evolve on timescales shorter than the structure persistence timescale associated 
with the long-lived flux tubes. 

Sheets and filaments are the dominant large-amplitude, long timescale structures that 
arise in the KAW system. Filaments arise and persist as long as rj is small, with their 
amplitude and statistical influence diminished as i] increases. Sheets exist in both regimes, 
becoming the sole large-scale structure in the large r] regime. Density gradients are consis- 
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tently non-Gaussian in both regimes as long as /i is small, although the density structures 
are different in both regimes. Density gradient sheets arise in the large t] regime and these 
density gradient sheets are large enough to yield non-Gaussian statistics. 

5. Ensemble Statistics and PDFs 

To explicitly analyze the extent to which the decaying KAW system develops non- 
Gaussian statistics, ensemble runs were performed for both the ////i ~ 1 and rj/fi 1 
regimes, and PDFs of the fields were generated. 

For the rj/fi ^ 1 regime, 10 numerical solutions were evolved with identical parameters 
but for different randomization seeds. In this case f] = fi and both damping parameters have 
minimal values to ensure numerical stability. The fields were initially phase-uncorrelated. 
The density gradient ensemble PDF for two times in the solution results is shown in Fig. 
Density gradients are Gaussian distributed initially. Many Alfven times into the numerical 
solution the statistics are non-Gaussian with long tails. These PDFs are consistent with 
the time histories of density gradient kurtosis excess as shown above. The distribution tail 
extends beyond 15 standard deviations, almost 90 orders of magnitude above a Gaussian 
best-fit distribution. Similar behavior is seen in the current PDFs - initially Gaussian 
distributed tending to strongly non-Gaussian statistics with long tails for later times. Fig. [T7] 
is the current PDF at an advanced time into the numerical solution. It is to be noted that the 
density gradient PDF has longer tails at higher amplitude than does the current PDF. One 
would expect these to be in rough agreement, since the underlying density and magnetic 
fields have comparable PDFs that remain Gaussian distributed throughout the numerical 
solution. The discrepancy between the density gradient and current PDFs suggests a process 
that enhances density derivatives above magnetic field derivatives. Future work is required 
to explore causes of this enhancement. This result is significant for pulsar scintillation, which 
is most sensitive to density gradients. Although interstellar turbulence is magnetic in nature, 
the KAW regime has the benefit of fluctuation equipartition between n and B. The density 
gradient, however, is more non-Gaussian than the magnetic component, suggesting that this 
type of turbulence is specially endowed to produce the type of scintillation scaling observed 
with pulsar signals. 

Ensemble runs for the rj/fi'^l regime yield distributions similar to the r^Z/i ~ 1 regime 
in all fields. The ensemble PDF for two times is shown in Fig. [181 The initial density gradient 
PDF is Gaussian distributed. For later times long tails are evident and consistent with the 
kurtosis excess measurements as presented above for the fi = case. The density gradient 
distribution has longer tails at higher amplitude than the current distribution; the overall 
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distributions are similar to those for the r^/yU ~ 1 regime, despite the absence of filamentary 
structures and the presence of sheets. The strongly non-Gaussian statistics are insensitive to 
the damping regime, provided that the diffusion coefficient is small enough to allow density 
gradients to persist. 

6. Discussion 

Using the normalizations for Eqs. ([3]) and (j4]) and using B = lAfiG, n = 0.08 cm~^ and 
Te = 1 eV, rjnorm, the normalized Spitzer resistivity, is 2.4 x 10~^ and finorm, the normalized 
coUisional diffusivity, is 1.9 x 10~^. For a resolution of 512^, these damping values are unable 
to keep the system numerically stable. The threshold for stability requires the simulation t] to 
be greater than 5 x 10~^, which is almost within an order of magnitude of the ISM value. The 
numerical solutions presented here, while motivated by the pulsar signal width scalings, more 
generally characterize the current and density gradient PDFs when the damping parameters 
are varied. We would expect the density gradients to be non-Gaussian when using parameters 
that correspond to the ISM. Future work will address the pulsar width scaling using electron 
density fields from the numerical solution. 

The non-Gaussian distributions presented here are strongly tied to the fact that the 
system is decaying and that circular intermittent structures are preserved from nonlinear 
interaction. One can show that, in the KAW system, circularly symmetric structures (or 
filaments) are force free in Eqs. ([3]) and (jlj), i.e., the nonlinearity is zero. Once a large- 
amplitude structure becomes sufficiently circularly symmetric and is able to preserve itself 
from background turbulence via the shear mechanism, that structure is expected to persist on 
long timescales relative to the turbulence. Structure mergers will lead to a time-asymptotic 
state with two oppositely-signed current structures and no turbulence. As structures merge, 
kurtosis excess increases until the system reaches a final two-filament state, which would 
have a strongly non-Gaussian distribution and large kurtosis excess. 

If the system were driven, energy input at large scales would replenish large-amplitude 
fluctuations. New structures would arise from large amplitude regions whenever the radial 
magnetic field shear were large enough to preserve the structure from interaction with tur- 
bulence. One could define a structure-replenishing rate from the driving terms that would 
depend on the energy injection rate and scale of injection. The non-Gaussian measures 
for a driven system would be characterized by a competition between the creation of new 
structures through the injection of energy at large scales and the annihilation of structures 
by mergers or by erosion from continuously replenished small-scale turbulence. If erosion 
effects dominate, the kurtosis excess is maintained at Gaussian values, diminishing the PDF 
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tails relative to a Levy distribution. If replenishing effects dominate, however, the enhance- 
ment of the tails of the density gradient PDF may be observed in a driven system as it is 
observed in the present decaying system. We note that structure function scaling in hy- 
drodynamic turbulence is consistent with the replenishing effects becoming more dominant 
relative to erosion effects as scales become smaller, i.e., the turbulence is more intermittent 
at smaller scales. The large range of scales in interstellar turbulence and the conversion of 
MHD fluctuations to kinetic Alfven fluctuations at small scales both support the notion that 
the structures of the decaying system are relevant to interstellar turbulen ce at the scales of 



KAW excitations. This scenario is consistent with arguments suggested by iHarmon fc Coles 



( 120051 ). They propose a turbulent cascade in the solar wind that injects energy into the KAW 
regime, counteracting Landau damping at scales near the ion Larmor radius. By doing so 
they can account for enhanced small-scale density fluctuations and observed scintillation 
effects in interplanetary scintillation. 

We also observe that, although the numerical solutions presented here are decaying in 
time, the decay rate decreases in absolute value for later times (Figs. [1] and |2]), approximating 
a steady-state configuration. The kurtosis excess (Figs. M and [TU]) for the density gradient 
field is statistically stationary after a brief startup period. Despite the decaying character of 
the numerical solutions, they suggest that the density gradient field would be non-Gaussian 
in the driven case. 

The kurtosis excess - a measure of a field's spatial intermittency - is itself intermittent 
in time. The large spikes in kurtosis excess correspond to rare events involving the merger 
of two large-amphtude structures, usually filaments. A large-amplitude short-hved sheet 
grows between the structures and persists throughout the merger, gaining amplitude in time 
until the point of merger. The kurtosis excess during this merger event is dominated by 
the single large-amplitude sheet between the merging structures. This would likely be the 
region of dominant scattering for scintillation, since a corresponding large-amplitude density 
gradient structure exists in this region as well. The temporal intermittency of kurtosis excess 
suggests that these mergers are rare and hence, of low probability. The heuristic picture of 
long undeviated Levy flights punctuated by large angular deviations could apply to these 
merger sheets. 



7. Conclusions 

Decaying kinetic Alfven wave turbulence is shown to yield non-Gaussian electron density 
gradients, consistent with non-Gaussian distributed density gradients inferred from pulsar 
width scaling with distance to source. With small resistivity, large-amplitude current fila- 
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ments form spontaneously from Gaussian initial conditions, and these filaments are spatially 
correlated with stable electron density structures. The electron density field, while Gaussian 
throughout the numerical solution, has gradients that are strongly non-Gaussian. Ensem- 
ble statistics for current and density gradient fields confirm the kurtosis measurements for 
individual runs. Density gradient statistics, when compared to current statistics, have more 
enhanced tails, even though both these fields are a single derivative away from electron den- 
sity and magnetic field, respectively, which are in equipartition and Gaussian distributed 
throughout the numerical solution. 

When all damping is placed in resistive diffusion (r//// — > regime), filamentary struc- 
tures give way to sheet-like structures in current, magnetic, electron density and density 
gradient fields. Kurtosis measurements remain similar to those for the small t] case, and the 
field PDFs also remain largely unchanged, despite the different large-amplitude structures 
at play. 

The kind of structures that emerge, whether filaments or sheets, is a function of the 
damping parameters. With rj and /i minimal to preserve numerical stability and of com- 
parable value, the decaying KAW system tends to form filamentary current structures with 
associated larger-scale magnetic and density structures, all generally circularly symmetric 
and long-lived. Each filament is associated with a flux tube and can be well separated from 
the surrounding turbulence. Sheets exist in this regime as well, and they are localized to the 
interface between flux tubes. With t] small and /i = 0, the system is in a sheet-dominated 
regime. Both regimes have density gradients that are non-Gaussian with large kurtosis. 

The effects on pulsar signal scintillation in each regime have yet to be ascertained 
directly. The conventional picture of a Levy flight is a random walk with step sizes distributed 
according to a long-tailed distribution with no defined variance. This gives rise to long, 
uninterrupted flights punctuated by large scattering events. This is in contrast to a normally- 
distributed random walk with relatively uniform step sizes and small scattering events. The 
intermittent filaments that arise in the small t] and fi regime are suggestive of structures that 
could scatter pulsar signals through large angles, however the associated density structures 
are broadened in comparison to the current filament and would not give rise to as large a 
scattering event. Ev en broadened structures can yield Levy distributed density gradients 



(ITerry fc SmithI 120071 ). but it is not clear how the Levy flight picture can be applied to these 
broad density gradient structures. In the /i = regime, the large-aspect-ratio sheets may 
serve to provide the necessary scatterings through refraction and may map well onto the 
Levy flight model. 

An alternative possibility, suggested by the temporal intermittency of the kurtosis (itself 
a measure of a field's spatial intermittency), is the encounter between the pulsar signal and 
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a short-lived sheet that arises during the merger of two filamentary structures. These sheets 
are limited in extent and have very large amplitudes. At their greatest magnitude they are 
the dominant structure in the numerical solution. Their temporal intermittency distinguish 
them from the long-lived sheets surrounding them. It is possible that a pulsar signal would 
undergo large scattering when interacting with a merger sheet. This scattering would be a 
rare event, suggestive of a scenario that would give rise to a Levy flight. 
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Fig. 1. — Energy vs. time for cross-correlated initial conditions. Total energy is monotoni- 
cally decreasing with time, and magnetic and internal energies remain in rough equipartition. 
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Fig. 2. — Energy vs. time for uncorrelated initial conditions. Total energy is monotonically 
decreasing with time, and magnetic and internal energies remain in rough equipartition. 
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Fig. 3. — Contours of n for various times in a numerical solution with correlated initial 
conditions. 




Fig. 4. — Contours of |B| for various times in a numerical solution with correlated initial 
conditions. 
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Fig. 5. — Contour plot of n with B vectors overlaid. The positive, circularly-symmetric 
density structures correspond to counterclockwise-directed B structures; the opposite holds 
for negative circularly-symmetric density structures. These spatial correlations are to be 
expected for correlated initial conditions. 
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Fig. 6. — Contours of n for various times in a numerical solution with uncorrelated initial 
conditions. 
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Fig. 7. — Contours of |B| for various times in a numerical solution with uncorrelated initial 
conditions. 
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Fig. 8. — Contour plot of n with B vectors overlaid for a numerical solution with initially 
uncorrelated initial conditions. The positive, circularly-symmetric density structures corre- 
spond to magnetic field structures, although the sense (clockwise or counterclockwise) of the 
magnetic field structure does not correlate with the sign of the density structures. Circled in 
black are symmetric structures that display a high degree of spatial correlation. The circle 
gives an approximate indication of the separatrix for the structure. 



-27- 




Fig. 9. — Kurtosis excess for a numerical solution with phase-correlated initial conditions 
and T]/ fj, = 1. 
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Fig. 10. — Kurtosis excess for a numerical solution with phase-uncorrelated initial conditions 
and T]/ fj, = 1. 
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Fig. 11. — Electron density contour visualization with diffusive damping parameter /i = 
for various times. 




Fig. 12. — Current density contour visualization with diffusive damping parameter jj, = 
for various times. 
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Fig. 13. — Magnitude of magnetic field contour visualization with diffusive damping param- 
eter fi = for various times. 
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Kurtosis excess vs. time -- zero-diffusivity 
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Fig. 14. — Kurtosis excess for a numerical solution with diffusive parameter /i = 0. Density 
gradient kurtosis remains greater than current kurtosis for the duration of the numerical 
solution. 
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Fig. 15. — Electron density gradient {x direction) contour visualization with diffusive damp- 
ing /i = for various times. 



- 33 - 



Log-PDF of Density Gradient witli Gaussian Best-Fit 




den. gradient, t=0 
den. grad. gaussian, t=0 
den. gradient, t=5000 
den. grad. gaussian, t=5000 



0.0 0.2 
Density gradient (dimensionless units) 



0.4 



0.6 



le- 



Fig. 16. — log-PDF of density gradients for an ensemble of numerical solutions with r]/ ^ = 1 
at if: = and t = 5000. The density gradient field at t = is Gaussian distributed, while for 
t = 5000 the gradients are enhanced in the tails, and deviate from a Gaussian. A best-fit 
Gaussian for each PDF is plotted for comparison. 
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Log-PDF of Current with Gaussian Best-Fit 




Fig. 17. — log-PDF of current for an ensemble of numerical solutions with r]/fi=latt = 
and t = 5000. The current at t = is Gaussian distributed. For t = 5000 the current is 
non- Gaussian. Unlike the density gradient, the current is not enhanced in the tails of the 
PDF for later times relative to its initial Gaussian envelope. 
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Fig. 18. — Log-PDF of density gradient for an ensemble of numerical solutions with /i = 
at if: = and t = 5000. The density gradient field at t = is Gaussian distributed, while for 
t = 5000 the gradients are enhanced in the tails, and deviate from a Gaussian. A best-fit 
Gaussian for each PDF is plotted for comparison. 



